# count average errors for each read between 150 and 320 bp
# and the number of reads with severe errors (>3%)

fa = ARGV[0]
xm = ARGV[1]

terr = 0

m3 = 0
m4 = 0

m5 = 0
m6 = 0
good = 1

ss = 0
total = 0
tlength = 0
File.new(fa,'r').each do |line|
  if line=~ /^>(\S+)/
    if  ss >= 150 and ss <= 320
      total += 1
      tlength += ss
    end
    ss = 0
  else
    ss += line.chomp.length
  end
end
if  ss >= 150 and ss <= 320
  total += 1
  tlength += ss
end

er = 0
ser = 0 

File.new(xm, "r").each do |line|
  if line=~ /^\s*(\d+)\s+(\d+\.\d+)\s+(\d+\.\d+)\s+(\d+\.\d+)\s+(\S+)\s+(\d+)\s+(\d+)\s+\((\d+)\)\s+/
    subs,ins,dels,name,qs,qe,qt = $2.to_f, $3.to_f, $4.to_f, $5, $6.to_i, $7.to_i, $8.to_i
    span = qe - qs + 1
    l = qe + qt
    er = ins + dels + subs
    next if l < 150 or l > 320
    terr += span*er
    if er >=6
      m6 += 1
    elsif er >=5 
      m5 += 1
    elsif er >= 4
      m4 += 1
    elsif er >=3
      m3 += 1
    elsif qt < 10 and qs < 10
      good += 1
    end
    if er >=3
      ser += span*er
    end
  end
end
ratio = terr / tlength
rp = ser / tlength
puts "#{fa}\t#{total}\t#{good}\t#{good/total.to_f}\t#{m3}\t#{m4}\t#{m5}\t#{m6}\t#{ratio}\t#{rp}"

    
